#!/bin/bash

b="-n 2 -m 1 --best --strata"
p=1
usage="
Author : Hyunmin Kim (hyun.kim@ucdenver.edu)
Tool   : Bowtie mapping in poly-seq style to the <out> directory 
Usage  : cat <fastq> | $0 [options] <bowtieindex> <out>
 <out>  : output directory (produces a.bam a.bam.bai)
 [options]   
    -b \"<str>\" : bowtie parameters (default: \"${b}\") 
    -p <int>  : number of processors (default: ${p}) 
    -h       : help
"
if [ -t 0 ] || [ $# -lt 1 ] || [ $1 == "-h" ]; then
    echo "$usage"
    exit 0
fi
while getopts "hp:b:" arg; do
    case $arg in
        p) p="${OPTARG}" ;;
        b) b="${OPTARG}" ;;
        h) echo "$usage"; exit 0 ;;
    esac
done
shift $(( OPTIND - 1 ));
if [ $# != 2 ]; then
    echo "Error: Not enough argument!!"
    exit 0
fi
bi=$1
out=${2%/}
mkdir -p $out

## run command
genome=${bi##*/}
cmd1="bowtie -p $p $b -S $bi - ";
cmd2=' chomp; my @a=split/\t/,$_; 
	if($a[0]=~/^@/){
		$a[1]=~s/dmel_mitochondrion_genome/M/;
		$a[1]=~s/SN:(\w+)/SN:chr$1/;
	}elsif($a[2] ne "\*"){
		$a[2]="chr".$a[2];
	}
	print join("\t",@a),"\n";
'

if [ ${genome} == "dm3" ]; then
	$cmd1 2> $out/log.bowtie | perl -ne "$cmd2"  | samtools view -bS - | samtools sort - $out/a
#| samtools view -bS - | samtools sort - $out/a
else
	$cmd1 2> $out/log.bowtie  | samtools view -bS - | samtools sort - $out/a
fi

#echo "$cmd";
#$cmd
#samtools index $out/a.bam
#"

#@SQ	SN:Uextra	LN:29004656
#@SQ	SN:2RHet	LN:3288761
#@SQ	SN:2LHet	LN:368872
#@SQ	SN:3LHet	LN:2555491
#@SQ	SN:3RHet	LN:2517507
